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[l]  One  difficulty  for  ocean  modeling  is  the  lack  of 
velocity  data  for  specifying  the  initial  condition.  Diagnostic 
initialization  is  widely  used;  it  integrates  the  model  from 
known  temperature  (Tc)  and  salinity  (Sc)  and  zero  velocity 
fields  while  holding  ( Tc ,  Sc)  unchanged.  After  a  period 
(around  30  days)  of  the  diagnostic  run,  the  velocity  field 
(Vc)  is  established,  and  (Tc,  Sc,  Vc)  fields  are  then  treated  as 
the  initial  conditions  for  the  prognostic  numerical  modeling. 
During  the  diagnostic  initialization  period,  the  heat  and  salt 
‘source/sink’  terms  are  generated  at  each  time  step. 
Maximum  time  rates  of  absolute  change  of  the  monthly 
mean  T,  S  (0. 1 7day,  0. 1  ppt/day)  are  taken  as  the  standard 
measures  to  identify  the  strength  of  the  thermohaline 
‘  sources/sinks ’.  Twenty  four  times  of  the  standard  measures 
(0.1°/hr,  0.1  ppt/hr)  represent  strong  ‘sources/sinks’.  Ten 
times  of  the  strong  ‘sources/sinks’  (l°/hr,  1  ppt/hr)  represent 
extremely  strong  ‘sources/sinks’.  The  Princeton  Ocean 
Model  implemented  for  the  Japan/East  Sea  is  used  to 
demonstrate  the  existence  of  extremely  strong  thermohaline 
sources  and  sinks  generated  by  the  diagnostic  initialization 
with  the  annual  mean  Tc,  Sc  from  the  Navy’s  Global  Digital 
Environmental  Model.  The  effects  of  extremely  strong  and 
spatially  nonuniform  initial  heating/cooling  (salting/ 
freshening)  rates  on  thermohaline  and  velocity  fields  need 
to  be  further  investigated.  INDEX  TERMS:  4263 
Oceanography:  General:  Ocean  prediction;  4255  Numerical 
modeling;  4243  Marginal  and  semienclosed  seas;  4520 
Oceanography:  Physical:  Eddies  and  mesoscale  processes;  4528 
Fronts  and  jets.  Citation:  Chu,  P.,  and  J.  Lan,  Extremely  strong 
thermohaline  sources/sinks  generated  by  diagnostic  initialization, 
Geophys.  Res.  Lett.,  30(6),  1341,  doi:10.1029/2002GL0 16525, 
2003. 

1.  Introduction 

[2]  Ocean  modeling  aims  to  integrate  hydrodynamic  and 
thermodynamic  equations  numerically  with  boundary  con¬ 
ditions  (lateral  and  vertical)  from  initial  states  of  temper¬ 
ature  (7),  salinity  ( S ),  and  velocity.  Initial  T,  S  fields  are 
relatively  easy  to  obtain,  such  as  using  climatological  ( Tc , 
Sc)  datasets  ([e.g.,  the  Navy’s  Global  Digital  Environmental 
Model  (GDEM)).  However,  the  initial  velocity  field  is 
usually  not  available  due  to  insufficient  number  of  velocity 
observations.  Thus,  initialization  of  the  velocity  field 
becomes  an  important  procedure  for  ocean  modeling. 
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[3]  A  widely  used  model  initialization  is  the  diagnostic 
mode,  which  integrates  the  model  from  known  T,  S  data 
such  as  climatological  ( Tc ,  Sc)  and  zero  velocity  fields, 
while  holding  ( Tc ,  Sc)  unchanged.  After  a  period  (about 
30  days)  of  the  diagnostic  run,  the  velocity  field  (Vc)  is 
established,  and  ( Tc ,  Sc,  Vc)  fields  are  treated  as  the  initial 
conditions  for  numerical  prognostic  modeling.  Since  initial 
condition  error  drastically  affects  model  predictability  [Lor¬ 
enz,  1963;  Chu,  1999],  questions  arise:  Does  the  diagnostic 
mode  provide  ideal  initialization?  What  is  the  physical 
process  associated  with  the  diagnostic  run? 

[4]  In  this  study,  the  Princeton  Ocean  Model  (POM) 
[Blumberg  and  Mellor,  1987]  implemented  for  the  Japan / 
East  Sea  (JES)  with  no  surface  and  lateral  forcing  is  used  to 
investigate  the  physical  outcome  of  the  diagnostic  initiali¬ 
zation.  The  GDEM  annual  mean  (Tc,  Sc)  data  with  0.5°  x 
0.5°  resolution  are  used. 


2.  Physical  Significance  of  the  Diagnostic  Mode 


[5]  Let  ( V,  w)  be  the  horizontal  and  vertical  velocity 
components,  and  V  the  horizontal  gradient  operator.  Ocean 
numerical  models  are  based  on  the  momentum  equation 


^  =  _V.W-»f-kx/V-iv,+  |(*r«f)  +  Hv. 


(1) 


and  the  temperature  and  salinity  equations, 

g.-v.vs-wg+lfyD+ft.  W 


where  /  is  the  Coriolis  parameter,  p  the  density,  p  the 
pressure,  and  (KM,  KH)  the  vertical  eddy  difiusivity  for 
turbulent  mixing  of  momentum,  temperature,  and  salinity. 
The  terms  (Hv,  Ht,  Hs)  represent  the  subgrid  processes 
causing  the  local  time  rate  of  change  in  (V,  T,  S). 

[6]  The  diagnostic  mode  procedure  integrates  (1)— (3) 
from 

T  =  Tc,  S  =  Sc,  V  =  0,  at  t  =  0,  (4) 

with  T  and  S  unchanged.  This  procedure  is  analogous  to  the 
process  of  adding  heat  and  salt  source/sink  terms  ( FT,  Fs)  in 
(2)  and  (3) 

dT  dT  d  („  dT\  „  „ 

—  =-V*VT-w— +  — (X/f— J  +Ht  +  Ft,  (5) 
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ds 

ft 


.V.VS-w|  +  |  (*«£)+*  +  *, 


(6) 


to  keep 


dT  „  dS  „ 
at  “  °’  dt  ~  ° 


(7) 


at  the  each  time  step.  Comparison  of  (7)  with  (5)  and  (6) 
leads  to 


fr,v.vr  +  »g-|(*r§)-«-.  (») 

F,  =  V.VS  +  .|-|(j6,|)-«,.  (9) 


Thus,  the  heat  and  salt  ‘source/sink’  terms  are  generated 
during  the  diagnostic  initialization  at  each  time  step. 


3.  Measures  of  ‘Source/Sink’  Strength 

[7]  Most  diagnostic  initialization  uses  the  climatologi- 
cally  monthly  (or  annual)  mean  data  as  the  initial  T,  S 
conditions.  The  maximum  variability  of  T,  S  is  estimated  by 
35°C  and  35  ppt.  Thus,  maximum  time  rates  of  absolute 
change  of  the  monthly  mean  T,  S  data  (usually  taken  as 
initial  conditions)  are  estimated  by 


dT 

dt 


K°f 

<  — -=*0.1°C  day-1, 

yr 


<  ~  0.1  ppt  day  *, 

yr 


(10) 

which  are  taken  as  the  standard  measures  for  ‘sources/ 
sinks’.  Twenty  four  times  of  the  standard  measures 


dT 

~0.1°C  hr-1, 

ds 

dt 

Strong 

dt 

~  O.lppt  hr  ’,  (11) 


Strong 


Figure  1.  Geography  and  isobaths  showing  the  bottom 
topography  of  the  JES  model. 

the  Okhotsk  Sea  through  the  Soya  and  Tatar  Straits.  The 
smoothed  JES  bathymetry  is  shown  in  Figure  1. 

[10]  The  Tsushima  Warm  Current  (TWC),  dominating  the 
surface  layer,  flows  in  from  the  Tsushima  Strait,  and  carries 
warm  water  from  the  south  up  to  40°N  where  a  subpolar 
front  (SPF)  forms.  Most  of  the  nearly  homogeneous  water 
in  the  deep  part  of  the  basin  is  called  the  Japan  Sea  Proper 
Water  [Moriyasu,  1972]  and  is  of  low  temperature  and  low 
salinity. 


represent  strong  ‘sources/sinks’.  Ten  times  of  the  strong 
‘sources/sinks’  \dT/dt\Snrmg  and  \dS/dt\strong, 


dT 

~  1°C  hr-1, 

ds 

dt 

Extra  Strong 

dt 

Extra  Strong 


~  lppt  hr  ’, 


(12) 


represent  extremely  strong  ‘sources/sinks’. 

4.  Diagnostic  Initialization  for  JES  Modeling 

[8]  Question  arises:  How  large  are  these  source/sink  terms 
after  the  initialization?  Are  the  false  ‘sources/sinks’  bearable 
for  numerical  modeling?  We  use  POM  implemented  for  the 
Japan/East  Sea  (JES)  to  evaluate  the  magnitude  of  the 
source/sink  terms  after  the  diagnostic  initialization. 

4.1.  Description  of  JES 

[9]  JES,  covering  an  area  of  106  km2,  has  steep  bottom 
topography  that  makes  it  a  unique  semi-enclosed  ocean 
basin  overlaid  by  a  pronounced  monsoon  surface  wind.  It 
has  a  maximum  depth  in  excess  of  3,700  m,  and  is  isolated 
from  open  oceans  except  for  small  (narrow  and  shallow) 
straits,  which  connect  the  JES  with  the  North  Pacific 
through  the  Tsushima/Korean  and  Tsugaru  Straits  and  with 


4.2.  Model  Description 

[11]  The  model  contains  94  x  100  horizontally  fixed  grid 
points.  The  horizontal  spacing  of  10'  latitude  and  longitude 
(approximately  11.54  to  15.18  km  in  the  zonal  direction  and 
18.53  km  in  the  latitudinal  direction)  and  15  vertical 
nonuniform  sigma  coordinate  levels.  The  model  domain  is 
from  35.0°N  to  52°N,  and  from  127.0°E  to  142.5°E.  The 
bottom  topography  is  obtained  from  the  smoothed  Naval 
Oceanographic  Office  Digital  Bathymetry  Data  Base  5 
minute  resolution.  The  horizontal  diflusivities  are  modeled 
using  the  Smagorinsky  form  with  the  coefficient  chosen  to 
be  0.2  for  this  application.  No  atmospheric  forcing  is 
applied  to  the  model. 

[12]  Closed  lateral  boundaries,  i.e.,  the  modeled  ocean 
bordered  by  land,  were  defined  using  a  free  slip  condition 
for  velocity  and  a  zero  gradient  condition  for  temperature  and 
salinity.  No  advective  or  diffusive  heat,  salt  or  velocity  fluxes 
occur  through  these  boundaries.  At  open  boundaries,  we  use 
the  radiative  boundary  condition  with  zero  volume  transport. 

4.3.  Mode  Splitting 

[13]  For  computational  efficiency,  the  mode  splitting 
technique  [Blumberg  and  Mellor,  1987]  is  applied  with  a 
barotropic  time  step  of  24  seconds,  based  on  the  Courant- 
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Figure  3.  Temporal  variation  of  total  kinetic  energy.  Note 
that  the  quasi-steady  state  is  reached  after  30  day’s 
integration. 

day-60  are  used  to  identity  their  magnitudes  and  sensitivity 
to  the  integration  period. 


Figure  2.  GDEM  annual  mean  (a)  temperature  and  (b) 
salinity  at  various  depths. 

Friederichs-Levy  (CFL)  computational  stability  condition 
and  the  external  wave  speed;  and  a  baroclinic  time  step  of 
720  seconds,  based  on  the  CFL  condition  and  the  internal 
wave  speed. 

4.4.  Climatological  T,  S  Data 

[m]  The  annual  mean  GDEM  T,  S  data  are  used  for  the 
study  [Chu  et  al.,  2001],  Here,  the  fields  at  the  surface  and 
150  m  depth  are  presented.  The  depth  of  150  m  corresponds 
to  the  permanent  thermocline  and  the  middle  level  of  the 
Japan  Sea  Intermediate  Water.  Below  the  depth  of  150  m,  the 
water  mass  is  quite  uniform.  Non-uniform  heat/salt  source 
and  sink  terms  might  cause  abrupt  thermohaline  change. 

[is]  The  climatological  mean  temperature  field  at  two 
depths  (0,  150  m)  clearly  shows  the  existence  of  the  SPF 
with  the  position  around  3  8°N  in  the  west  and  near  42°N  in  the 
east  (Figures  2a  and  2b).  The  temperature  is  more  than  6°C 
higher  south  of  SPF  than  north  of  SPF  at  0  and  1 50  m.  The 
climatological  mean  salinity  field  at  the  surface  (Figure  2c) 
clearly  shows  that  the  saline  Kuroshio  water  (>34.2  ppt) 
enters  the  JES  from  the  Tsushima/Korean  Strait  into  the  JES 
and  forms  two  permanent  salty  centers  with  the  salinity  higher 
than  34.0  ppt,  located  north  of  SPF  (west  of  the  Hokkaido 
Island)  and  south  of  SPF  at  37°-40°N,  132°-136°E. 

4.5.  Diagnostic  Mode 

[i6]  The  POM  was  integrated  in  the  diagnostic  mode 
with  all  three  components  of  velocity  («,  v,  w)  initially  set  to 
zero,  and  with  temperature  and  salinity  specified  by 
interpolating  GDEM  annual  mean  data  to  each  model  grid 
point.  Ft  and  Fs  are  obtained  at  each  time  step.  The 
diagnostic  model  was  integrated  for  60  days,  30  days  were 
sufficient  for  the  volume-mean  model  kinetic  energy  to 
reach  quasi-steady  state  under  the  imposed  conditions 
(Figure  3).  The  thermohaline  source/sink  terms  (Fr,  Fs) 
generated  by  the  diagnostic  initialization  on  day-30  and 


5.  Source/Sink  Terms  in  Diagnostic  Initialization 

5.1.  Heat  Sources/Sinks 

[17]  Horizontal  distributions  of  pCpFT  (unit:  W  m-3)  on 
day-30  (Figure  4)  at  the  four  <7  levels  (0,  —0.143,  —0.5,  and 
—0.929)  show  extremely  strong  heat  sources/sinks  gener¬ 
ated  by  the  diagnostic  initialization.  For  pc^r  =  1000  W 
m-3,  the  time  rate  of  absolute  temperature  change  Fr  = 
0.84°C  hr  "1.  The  sources/sinks  have  various  scales  and 
strengths.  They  reveal  small-  to  meso-scale  patterns  in  most 
areas  except  a  large-scale  pattern  near  the  bottom  (ct  = 
-0.929).  The  strength  of  the  source/sink  increases  with 
depth  from  the  surface  to  subsurface.  The  extremely  strong 
source  reaching  5 1 64W  m~3  (corresponding  to  Fr  =  4.34°C 


Longitude(E) 


Longitude(E) 


Figure  4.  Horizontal  distribution  of  pc^F,  (unit:  W  m  3) 
on  day-30  at  a  levels  of:  (a)  0,  (b)  -0.143,  (c)  -0.5,  and  (d) 
-0.929. 
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Figure  5.  Horizontal  distribution  of  Fs  (unit:  ppt  s  ')  on 
day-30  at  a  levels  of:  (a)  0,  (b)  -0.143,  (c)  -0.5,  and  (d) 
-0.929. 

hr  ')  and  the  extremely  strong  sink  reaching  —4735  Win  3 
(corresponding  to  Fr=  — 3.98°C  hr-1),  and  decreases  with 
depth  below  the  subsurface.  Near  the  bottom,  the  JES  basin 
is  dominated  by  cooling  with  the  maximum  sink  strength 
—  1557  Wm-3  (corresponding  to  Fr=  —  1.31°C  hr-1).  From 
the  subsurface  to  the  bottom,  the  source/sink  terms  have 
some  organized  pattern  near  the  SPF.  At  the  subsurface,  a 
dipole  pattern  occurs  between  133°-136°E  with  strong 
source  strength  approximate  2000  W  m-3  (corresponding  to 
Ft  =  1.68°C  hr-1)  north  of  the  SPF  and  cooling  rate  (~ 
—2000  W  m-3)  south  of  the  SPF.  Near  the  bottom,  a  large 
cooling  area  with  the  cooling  rate  of  750  Wm-3  occurs 
north  of  the  SPF. 

5.2.  Salt  Sources/Sinks 

[is]  Horizontal  distributions  of  Fs  (unit:  ppt  m-3)  on  day- 
30  (Figure  5)  at  the  four  a  levels  (0,  —0.143,  —0.5,  and 
—0.929)  show  near-extremely  strong  salinity  sources/sinks 
generated  by  the  diagnostic  initialization.  These  sources/ 
sinks  have  various  scales  and  strengths.  They  reveal  small- 
to  meso-scale  patterns  in  most  areas  but  a  large-scale  pattern 
in  the  southern  JES  near  Tsushima/Korean  Strait  at  the 
surface  and  north  of  the  SPF  at  the  mid-level  and  bottom. 
The  strength  of  the  source/sink  increases  with  depth  from 
the  surface  to  the  bottom.  The  maximum  salinity  source 
(sink)  is  found  0.80  ppt  hr-1  (—0.76  ppt  hr-1)  at  cr  = 
-0.929  (a  =  -0.5). 

6.  Discussion 

[19]  When  the  prognostic  integration  starts,  the  source/ 
sink  terms  FT  and  Fs  are  removed  from  (5)  and  (6). 
Extremely  strong  and  spatially  nonuniform  initial  heating/ 
cooling  (salting/freshening)  rates  are  introduced  in  the 


ocean  models  and  cause  drastic  changes  in  thermohaline 
and  velocity  fields  initially  (after  the  diagnostic  run) 
especially  in  the  deep  layer  below  the  thermocline  and 
halocline.  Note  that  the  problem  is  caused  by  the  diagnostic 
initialization  only,  nothing  to  do  with  the  ocean  model  itself. 

[20]  In  the  diagnostic  initialization,  the  source/sink  terms 
drive  the  velocity  through  the  pressure  gradient  force  (see 
equation  (1)).  The  pressure  gradient  error  leads  to  errors  in 
the  initialized  velocity  field.  Different  models  (z-level,  ct- 
level,  and  5-level)  have  different  pressure  gradient  errors, 
which  in  turns  generate  different  initial  velocity  fields. 
Besides,  the  diagnostic  process  (spin  up/down)  largely 
depends  on  the  diffusion.  The  spin  down  scale  of  30  days  is 
the  state  of  balance  between  the  pressure  gradient  force  (not 
change  with  time)  and  other  terms  in  equation  (1)  such  as 
the  diffusion  term  that  depends  not  only  on  the  velocity 
field,  but  also  on  the  model  parameters.  Thus,  the  diagnostic 
initialization  depends  on  model  type  and  model  parameters. 

[21]  If  the  diagnostic  initialization  continues  to  be  used,  it 
is  urgent  to  study  the  following  problems:  Does  this 
artificial  initial  heating/cooling  (salting/freshening)  induce 
false  chaotic  motion  in  ocean  models?  How  long  does  the 
ocean  model  need  to  be  adjusted?  Does  the  spin-up  of  the 
prognostic  mn  have  the  capability  to  diminish  this  initial 
effect? 

[22]  If  the  monthly  mean  Tc,  Sc  data  are  used  as  the  initial 
conditions,  the  initial  heating  (or  cooling)  and  salting  (or 
freshening)  rates  should  not  be  greater  than  the  standard 
measures  (10)  everywhere  in  the  domain.  If  they  reach  the 
levels  of  strong  ‘sources/sinks’  ( \dTldt\Strong  and  |95/ 
dAstrong),  the  calculated  ( Tc ,  Sc,  Vc)  fields  after  diagnostic 
initialization  are  abnormal.  If  they  reach  the  levels  of  extra 
strong  ‘sources/sinks’  {\dTldt\ExtraStrong  and  | dS/ 
dt  I  ExtraStron g)i  the  calculated  ( Tc ,  Sc,  Vc)  fields  cannot  be 
used.  Thus,  development  of  a  check-up  algorithm  on 
strength  of  the  initial  source  and  sink  is  urgent. 
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